Correlation effects in the ground state charge density of Mott-insulating NiO: a 
comparison of ab-initio calculations and high-energy electron diffraction measurements 
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Accurate high-energy electron diffraction measurements of structure factors of NiO have been 
carried out to investigate how strong correlations in the Ni 3d shell affect electron charge density in 
the interior area of nickel ions and whether the new ab-initio approaches to the electronic structure 
of strongly correlated metal oxides are in accord with experimental observations. The generalized 
gradient approximation (GGA) and the local spin density approximation corrected by the Hubbard 
U term (LSDA+U) are found to provide the closest match to experimental measurements. The 
comparison of calculated and observed electron charge densities shows that correlations in the Ni 
3d shell suppress covalent bonding between the oxygen and nickel sublattices. 

PACS No: 61.14.Lj, 71.15.Mb, 71.20.Be 

Recent years have witnessed the largely unexpected |lj] 
progress in the development of computational approaches 
to the evaluation of fundamental properties of materials 
from the first principles. The stimulus for this develop- 
ment was provided by the Hohenberg-Kohn theorem ||] 
which establishes that the energy of the ground state of 
a solid is a functional of its one-electron density p(r). 
The problem of accurate determination of p(r) therefore 
acquires fundamental significance for the physics of mate- 
rials both from the experimental and theoretical points of 
view. In cases where accurate experimentally measured 
and calculated charge densities are available (like, e.g. in 
the case of silicon M), the quality of ab-initio approxi- 
mations can be assessed on the basis of the agreement 
between experimental and theoretical data. 

The Kohn-Sham method M, which provides a conve- 
nient way of carrying out density functional calculations, 
in certain cases encounters serious difficulties. For ex- 
ample, it predicts metallic ground states for a number of 
late transition metal monoxides where metal ions have 
partly filled electronic shells. Nickel and cobalt monox- 
ides are often quoted H as typical examples illustrating 
the failure of conventional density functional methods 
to describe the effective one-particle band structure of 
Mott insulating materials [glJ^] ■ Several modified density 
functional schemes have been proposed lately to explain 
the nature of large bandgaps observed for CoO and NiO. 
These schemes include the orbital polarization correction 
pL the self-interaction correction (SIC) M and the lo- 
cal spin density approximation including the Hubbard U 
term (LSDA+U) |H],[Tl|]. The new approximations im- 
prove the description of the effective one-particle band 
structure of Mott insulators, and their validity is fur- 
ther confirmed by the recent studies of orbital ordering 
in transition metal compounds, see e.g. ]l2| |. 

At the same time it is widely appreciated that the new 



'improved' computational schemes represent a departure 
from the original formulation of density functional theory 
g . The new approximations employ functionals that de- 
pend not only on the spin density of electrons p a (r) but 
also on the orbital occupation numbers that in turn de- 
pend on the choice of the basis functions. Calculations 
performed using the new ab-initio schemes result in bet- 
ter values for bandgaps and magnetic moments pQ ], At 
the same time the use of the modified energy functionals 
alters the relative occupancies of d-states [|l3| and alters 
the distribution of the charge density in the unit cell. 
This raises the question of how well the new function- 
als describe the main entity of density functional theory, 
namely, p(r) itself. 

In this paper we investigate this issue by comparing the 
calculated and experimentally observed charge density 
distributions. We compare the structure factors of NiO 
that were measured using a recently developed highly 
accurate electron diffraction technique |T3-n(| and cal- 
culated theoretically using several ab-initio linear muffin- 
tin orbital (LMTO)-based methods 17 lq], including the 
LSDA+U approach. 

At present, there is no sufficiently accurate experimen- 
tal information on the distribution of electron charge 
density in the unit cell of NiO or other similar transi- 
tion metal oxides. The powder X-ray diffraction tech- 
niques that were successfully used to determine the equi- 
librium positions of ions in a unit cell |19Jj20(| , do not have 
the accuracy required for observing the relatively small 
changes in the charge density resulting from the compe- 
tition between covalent bonding and correlation effects. 
The convergent beam electron diffraction (CBED) tech- 
nique used here takes advantage of the fact that electron 
beam can be focused on a small nearly perfect area of 
the sample and the resulting diffraction pattern can be 
simulated using highly accurate multiple scattering dy- 



namical diffraction approach |T|] therefore eliminating 
the extinction problem that limits the accuracy of X- 
ray techniques. The high precision of electron diffraction 
measurements has made it possible to study subtle de- 
tails of the charge density distribution in band insulators 
like MgO @ and Cu 2 Q. In this paper we for the 
first time investigate how the interplay between covalent 
bonding and the Coulomb on-site repulsion between d- 
electrons in a partially filled shell influences the charge 
density distribution in a crystal unit cell of a Mott insu- 
lating oxide. 

Nickel monoxide is probably the most extensively stud- 
ied Mott insulating material |2^,g3|] and it is often re- 
ferred to as the prototype of the entire class of 'anoma- 
lous' transition metal oxides. Depending on the type 
of approximation used in an ab-initio calculation, the 
ground state of NiO is predicted to be a metal (LDA, 
non-magnetic state), a ^0.4 eV bandgap Mott insula- 
tor (LSDA, antiferromagnetic state), a ~1.0 eV bandgap 
Mott insulator (GGA, antiferromagnetic state) or a ~3.0 
eV bandgap charge-transfer insulator (LSDA+U, antifer- 
romagnetic state). 



Figure 1 near here 



Fig.l shows plots of the density of states calculated 
for the above four cases. The X-ray photoelectron spec- 
troscopy data [£3| agree best with the LSDA+U one- 
particle band structure |H],g4| that shows that NiO 
is a charge-transfer insulator where the bandgap sepa- 
rates filled oxygen 2p and empty nickel 3d states. The 
band structures of NiO calculated using either LSDA or 
GGA show instead that the bandgap separates filled and 
empty nickel 3d states and that NiO is therefore a Mott- 
Hubbard insulator. 

To characterize the distribution of electron charge den- 
sity in a unit cell of NiO, we measured seven low order 
energy dependent structure factors [£5| U(G) that are 
defined by 
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where G is a reciprocal lattice vector, m is the relativistic 
electron mass and f2 is the volume of the unit cell. Sum- 



mation over a is performed over ions in a unit cell, Z a is 
the charge, r a is the equilibrium position and u a is the 
thermal displacement of the respective nucleus. p( T \r) is 
the electron density averaged over the thermal ensemble. 

The experiment was performed using the LEO-912 fi 
energy-filtering electron microscope with the Gatan liq- 
uid nitrogen cooled sample holder. The specimen used 
is a single crystal NiO cooled to about 110K. The small 
rhombohcdral distortion of NiO at 110 K was measured 
using higher-order Laue zone lines fl26| , f27[ to be a — 4.25A 
and a = 90.044°. This small distortion was neglected 
in the charge density study. The experimental CBED 
patterns were recorded using a 15 eV energy-filtering 
slit that was placed around the zero-loss peak. This 
was done to remove the contribution from inelastically 
scattered electrons that form continuous background due 
to plasmon and higher energy loss processes. Off-zone- 
axis systematic diffraction conditions were used to collect 
diffraction intensities for reflections up to (440). The ex- 
perimental patterns recorded using a slow-scan charge 
coupled device (CCD) camera were processed for the sub- 
sequent fitting using procedures described in p5[ . The 
refined values of the structure factors were obtained us- 
ing the EXTAL program |ljj. Fig. 2 illustrates the level 
of agreement between the measured intensity variations 
and multiple scattering dynamical diffraction simulations 
used in the refinement procedure. The error in the mea- 
sured structure factors was estimated by comparing the 
results obtained using line scans taken at different posi- 
tions (see Fig. 2). 

The calculated values of structure factors (0) were 
obtained by integrating the self-consistent solutions 
p(r) of the Kohn-Sham equations found using the 
LMTO method p7| and various approximations for the 
exchange-correlation potential. Calculations were per- 
formed using three k— panels and 512 k-points in the 
Brillouin zone. Convergence of the calculated values of 
structure factors was ensured by varying the number of 
k-points and by introducing additional approximations 
(e.g. by taking into account the spin-orbit coupling). 
The exchange-correlation functionals were taken from 
@ (LSDA), || (GGA) and || (LSDA+U). The full 
self-consistent charge density was represented by a sum of 
two terms [fl8| where the first term p(r) approximates the 
density in the region between the muffin-tin spheres and 
is continuous across the boundaries of the spheres. The 
second term approximates the density inside the muffin- 
tine spheres and is represented by a spherical harmonics 
expansion Yi m (9, </>). Substituting this in p) we arrive at 
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where ji(Gr) is a Bessel function, and u denotes the 
effective amplitude of thermal vibrations characterizing 
the motion of electrons in the interstitial region. No 
other thermal effects are taken into account in (g) in ac- 
cord with spectroscopical data J31J showing no detectable 
changes in the electronic structure of NiO occurring in 
the temperature range between 0°K and 615°K. 

A meaningful comparison between experiment and the- 
ory has to take into account the thermal vibrations of 
atomic nuclei. For example, the value of the (111) 
structure factor calculated in the LSDA(AF) approxi- 
mation assuming that nuclei are frozen in their equilib- 
rium positions equals —4.488 • 10~ 2 A 2 to be compared 
with —4.597 ■ 10 -2 A 2 that was obtained assuming that 
B Nl = 0.131 A 2 and B = 0.239 A 2 . The difference be- 
tween the two values is approximately ten times the ex- 
perimental uncertainty in the determination of this struc- 
ture factor (see Table) and this illustrates the significance 
of taking thermal vibrations into account. However, no 
reliable independent X-ray or neutron diffraction mea- 
surements on the Debye- Waller factors of Ni and O ions 
is available in the literature (notably, negative values of 
Bni were reported in a recent publication p0fl ). To pro- 
vide a starting approximation for a subsequent refined 
search, we calculated the temperature factors of Ni and 
O using the shell model. 

The shell model used for the theoretical evaluation of 
the Debye- Waller factors was developed by the Chalk 
River group [p2[ . In this model an ion is represented 
by a massive core and a rigid shell describing valence 
electrons. Nine parameters in total were introduced to 
describe the lattice dynamics of NiO, these include the 
shell charges, force constants for springs connecting the 
cores and the shells as well as the first nearest neigh- 
bours (for the Ni 2+ ions) and up to the second nearest 
neighbours (for the 2 ~ ions). The model parameters 
were obtained by fitting the calculated phonon dispersion 
curves to the experimentally measured ones p3| . For the 
chosen set of parameters of the model, we calculated the 
average thermal displacements (u 2 ) for all the modes of 
lattice vibrations, and also the Debye- Waller factors for 
both the nickel and oxygen ions. Fig. 3 shows the fitted 
phonon dispersion curves plotted for the (111) and (110) 
directions, and also the temperature dependent Debye- 
Waller factors. We have also investigated several other 
implementations of the shell model but found that they 
led to no significant improvement in the description of 
the phonon dispersion curves [p4| . 

To find more accurate values of the Debye- Waller pa- 
rameters for each choice of the exchange-correlation po- 



tential used for ab-initio calculations we plotted two- 
dimensional maps of the reliability factor R treating Bj^i 
and Bq as independent variables. 



Figure 4 near here 



The values of Bn% and Bo corresponding to the min- 
imum of the i?-factor were then used to obtain the val- 
ues of U (G) shown in the Table. Results listed in the 
Table show that the estimated values of Bjsh and Bo 
are nearly independent on the choice of approximation 
used in ab-initio calculations and that the spread of val- 
ues of the Debye- Waller factor does not exceed 6%. The 
value of the Debye- Waller factor characterizing the ther- 
mal motion of electrons in the interstitial region was eval- 
uated using two different approximations, namely, (B) = 
(B +B Ni )/2oi (B) = (M Bo+M m B m )/(M +M m ). 
The difference between structure factors evaluated using 
these two approximations was found to be significantly 
smaller than the uncertainty of experimentally measured 
values of structure factors. 

Apart from values calculated using the superposition 
of atomic densities, all the ab-initio methods exhibit high 
(better than 1%) degree of accord with experimental 
data, with the exception of (111) structure factor. There 
is a large spread among theoretical values of the (111) 
structure factor, which increases significantly with the 
inclusion of U. The structure factor of (111) is most sen- 
sitive to the changes in the distribution of the density of 
valence electrons, and the differences among the theoret- 
ical models shows primarily the differences in the calcu- 
lated ground state valence charge density. In terms of 
the overall R-factor, the closest approximation to the ex- 
periment is provided by the generalized gradient approx- 
imation and by the LSDA+U approach (see Table). Bet- 
ter agreement with the experimental value of the (111) 
structure factor is achieved with the GGA. 

The difference between the GGA and the LSDA+U 
charge density distributions is illustrated in Fig. 5 where 
we mapped the densities calculated using the GGA and 
the LSDA+U approximations, subtracting from each of 
them the density corresponding to the non-magnetic 
LDA solution. 



Figure 5 near here 



Fig. 5 shows that the symmetry of the deformation of 
electron density resulting from correlation effects remains 
the same both in LSDA+U and in GGA. At the same 
time there are significant differences in the radial struc- 
ture of the density distributions around Ni ions calcu- 
lated for the two cases. The LSDA+U approximation 
treats the wave functions of 3d states as 'rigid' objects 
where the Hubbard correction shifts the filled and unoc- 
cupied states in the opposite directions along the energy 
axis (2jj3^]. In the GGA approximation the shape of 
wave functions and the population of the 3d states de- 
pends on the local density and its gradient in the interior 
area of nickel ions. The LSDA+U approximation relies 
to a larger extent on the model assumptions and on the 
choice of tight-binding orbitals used for treating electron 
correlations in a partly filled 3d shell. The GGA approxi- 
mation uses the one-electron orbitals as auxiliary entities 
required in a calculation of total density p(r), which is 
the quantity observed experimentally using high-energy 
electron diffraction. 

Fig. 5 also shows a low resolution difference map be- 
tween the experimentally observed and the GGA charge 
density distribution estimated using seventy six low-order 
experimentally measured structure factors listed in the 
Table (this includes transpositions and mirror reflec- 
tions). The comparison between the experimental and 
calculated distributions confirms the trends revealed by 
the GGA and LSDA+U analysis showing that correla- 
tion effects arc responsible for the suppression of cova- 
lent bonding between the metal and oxygen sublattices 
(this effect manifests itself in the reduction of the charge 
density in the areas between the oxygen and nickel ions). 
This agrees with the analysis of a similar effect discovered 
in |3q| for uranium dioxide. The charge density in a unit 
cell of real NiO is more concentrated around atomic nu- 
clei and it also shows tendency towards increasing in the 
region between ions of the same type. Qualitatively this 
may be interpreted as an indication that Ni rf-orbitals in 
fact have the shape that is different from that predicted 
by either the GGA or the LSDA+U calculations. 

There are several reasons responsible for the observed 
inaccuracy of ab-initio approaches. For example, the 
LSDA+U approximation is based on the mean-field treat- 
ment of correlation effects [M. There could be other, 
more fundamental, circumstances leading to the differ- 
ence between the calculated and experimentally observed 
charge density distributions. One of the ideas behind the 
development of more accurate approaches to the treat- 
ment of electron correlations in transition metal com- 
pounds consists in that the new approaches are intended 



to be used for evaluating the parameters entering tight- 
binding many-body models of electron-electron interac- 
tions. These models are always based on a particular 
choice of orbitals associated with each of the ions in the 
solid. Our results show that the accuracy of the assump- 
tion that the charge density may be decomposed into 
contributions associated with individual ions, is limited, 
and this conclusion agrees with the analysis performed in 
Rcf. [Q. To describe correlation effects in oxides where 
covalency as well as correlation effects play a significant 
part, it may be necessary to take into account changes 
in both the shape and occupation of localized electronic 
orbitals. A second-quantized many-body model describ- 
ing intersite hopping and on-site Coulomb interaction be- 
tween electrons may prove to be sufficient for accounting 
for the positions of the main peaks in the spectrum of 
excited states of an oxide. At the same time even the 
exact solution of the model may not be capable of giv- 
ing a sufficiently accurate description to the the ground- 
state properties of the system such as the distribution 
of the charge density of electrons in a unit cell. The ap- 
proach developed in this paper can also now be applied to 
test the accuracy and to compare several other ab-initio 
methods that we did not consider above, for example, 
the self-interaction correctio n |9 ] or the first-principles 
Hartree-Fock approximation |36[ |. 

In summary, by combining a recently developed high 
accuracy electron diffraction technique with ab-initio cal- 
culations, we investigated how electron correlations in 
the Ni 3d shell influence the distribution of charge den- 
sity in the unit cell of NiO. By comparing the experi- 
mentally measured values of structure factors with values 
calculated using several different ab-initio approaches we 
found that the structure factors evaluated using the gen- 
eralized gradient approximation and the LSDA+U ap- 
proach agree best with the available experimental infor- 
mation. The experimental data show that the degree of 
covalent bonding in NiO is smaller than that predicted 
by theoretical calculations. 
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FIGURE CAPTIONS 



FIG. 1. The total density of states plots of NiO calculated 
assuming (a) a non-magnetic ground state and the local den- 
sity approximation (b) the type II AF ground state and the 
local spin density approximation (c) the type II AF ground 
state and the generalized gradient approximation (d) the type 
II AF ground state and the local spin density approximation 
corrected by the Hubbard U term (LSDA+U). ef denotes 
the Fermi energy and e e ^ p shows the experimentally observed 
position of the bottom of the conduction band. 



FIG. 2. An example of NiO structural factor measurement 
using convergent beam electron diffraction (CBED). The top 
is the experimentally recorded diffraction pattern with (220) 
and (440) strongly diffracted. The structural factors of (220) 
and (440) were obtained by fitting intensities along the in- 
dicated lines. The best fit is shown in bottom right. The 
schematic diagram on bottom-left shows the formation of 
CBED by focusing the electron beam on the top of crystalline 
specimen. 



FIG. 3. Plots illustrating the fitting of the phonon disper- 
sion curves along (111) and (1 10) directions, and the compar- 
ison between the values of the Debye- Waller factors evaluated 
using the shell model and found by comparing the experimen- 
tally measured and calculated structure factors. 



FIG. 4. A map showing the dependence of the R-factor on 
the two Debye- Waller factors B^i and Bo characterizing the 
amplitude of thermal vibrations of Ni and O ions in NiO. 
Values B Ni - 0.133 A 2 and B Ni - 0.244 A 2 correspond to 
the minimum of R = 0.0049. 



FIG. 5. Cross-sections of charge density distribution in the 
(100) plane of NiO calculated using the GGA and LSDA+U 
approximations. The two contour maps on the left-hand side 
show the difference between the self-consistent density distri- 
butions and the density calculated using the local density ap- 
proximation for a non-magnetic ground state. The map on the 
right-hand side shows the difference between the experimen- 
tally observed charge density distribution and the distribution 
calculated using the generalized gradient approximation. 



Table. The observed and calculated values of structure factors for NiO. The Debye- Waller temperature factors 
Bni = 87r 2 (u^ ri ) and Bo = 8ir 2 {u ) for each calculated set were introduced following the procedure described in the 
text. The energy of the incident electrons equals Eq = 119.5 keV. All the values listed in the table are given with the 
opposite sign and in A 2 units (see Ref. [||). The i?-factor is defined as R = £ G W G \U th (G) - [/ oxp (G)|/|£/ cxp (G)|, 
where U th (G) are the calculated and C/ cxp (G) are the experimentally measured values. The weight factors Wg are 
given by Wg = <7q /(X)g °g ) wnere "G represent experimental uncertainties. The error-bar of the i?-factor, 6R, is 
given by (SR) 2 — J2g^ g<t g/\^ cxp (^')\ 2 - Abbreviations NM and AF refer to non-magnetic and antiferromagnetic 
states, respectively. 
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-[/(111) 


4.632 


10" 2 (±0.012 


io-") 


4.401 


10 -2 


4.555 


io-" 


4.597 


io-" 


4.668 


io-" 


4.642 


10 -2 


4.708 


io-" 


-t/(200) 


9.083 


10"" (±0.022 


10"") 


9.486 


10" 2 


9.187 


io-" 


9.204 


io-" 


9.181 


10"" 


9.173 


10" 2 


9.151 


10"" 


-U(220) 


6.640 


10 -2 (±0.036 


io-") 


6.756 


IO"" 1 


6.709 


io-" 


6.716 


io-" 


6.705 


io-" 


6.703 


IO"" 1 


6.694 


10"" 


-Z7(311) 


2.482 


10 -2 (±0.010 


io-") 


2.482 


io- y 


2.482 


io-" 


2.482 


io-" 


2.481 


io-" 


2.481 


io- y 


2.481 


io-" 


-U(222) 


5.187 


10 -2 (±0.026 


10"") 


5.336 


10" 2 


5.325 


io-" 


5.318 


io-" 


5.307 


10"" 


5.311 


10" 2 


5.302 


10"" 


-U(4Q0) 


4.456 


10 -2 (±0.018 


10"") 


4.427 


10" 2 


4.455 


io-" 


4.456 


io-" 


4.457 


10"" 


4.456 


10" 2 


4.457 


10"" 


-C/(440) 


2.614 


10 -2 (±0.034 


10-") 


2.613 


10" y 


2.613 


IO - " 


2.611 


IO - " 


2.603 


10"" 


2.610 


10" y 


2.603 


10"" 


R, SR 


0.0045**, 0.0018 


0.0215 


0.0085 


0.0067 


0.0064 


0.0049 


0.0078 



values calculated using the shell model (see text) 
* R-factor evaluated on the basis of experimental uncertainties ug 
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I L»l I I I I I I | I I I I I I I I I 

1 LDA(NM)) 



* i i i i i i i i | i i i i i i i i 

LSDA(AF) 



i i i i i i i 



i i i i i 



GGA(AF) 



i i i i i i i i i i i i i i i i i i 

LSDA+UE 




.exp 



density of states (number of electrons per unit cell per eV) 



0.35 



q=2rc(i;,£,0)/a 



q=2it(ttQ/a 



0.10 



O 0.30 


— 






** 




(A 




O) 




C 




< 




w ° 25 




l_ 




o 




■1— » 




o 




CO 




M— 




<5 0.20 




«J 




Sf 




a> 




>• 0.15 




a> 




a 


.-*- " 



A GGA-Ni 

▼ GGA-0 
Shell model - Ni 
Shell model - O 



65 



115 165 215 

Temperature (K) 



a Experimental 
Shell model 




LO 


^^-^To 






-^- TO 




/"" i ~^\LA 


LA , 






X±^jX 






aS, \ 




-^TA 


T 2 A N. \ 

1 , 1 A 


f 1 


1 , 1 



20.0 



16.0 



12.0 2 

C 

o 

3 
O 

«< 



- 8.0 



- 4.0 



N 



265 -0.! 



-0.4 



-0.2 0.0 0.2 0.3 0.5 



0.0 




0.11 0.12 0.13 0.14 0.15 0.16 



%i 



ABOVE 
0.008 - 
0.006 - 
0.004 - 
0.002 - 
0.000 - 
-0.002 - 
-0.004 - 
-0.006 - 
-0.008 - 
-0.010- 
BELOW 



0.010 

0.010 

0.008 0.9 

0.006 

0.004 

0.002 

0.000 

-0.002 0.7 
-0.004 
-0.006 6 
-0.008 
-0.010 




o.o-p 



0.0 0.1 0.2 0.3 0.4 0.5 0.0 0.1 0.2 0.3 0.4 0.5 0.0 0.1 0.2 0.3 0.4 0.5 



